Analysis of an Eucalyptus clonal plantation dataset
In this analysis we test whether intraspecific variability (IV) in observed individual tree growth can emerge from the environment. To this aim, we use a clonal experimental setup. The EUCFLUX common garden is a clonal trial where 14 Eucalyptus clones are grown in a replicated, statistically-sound design. One of its main goals is to determine and compare the productivity of each clone. Our hypothesis on growth IV is that it mainly results from responses to environmental factors, and not only from intrinsic (here more precisely genetic) factors. Therefore, we use the EUCFLUX dataset in order to quantify IV of growth within clones, i.e. in a dataset where genetically-driven IV is nil. Following our hypothesis, we expect to detect IV in growth within single genotypes, although we also expect the genotype to influence the growth response.
1 The dataset
The EUCFLUX experiment is located in Brazil, in the state of São Paulo. It includes 14 clones of 5 different species or hybrids, planted at the same date and grown in 10 replicated blocks of 100 trees each, which were monitored during 6 years. Each replicated block is distant by up to 1.5 km in a 200 ha stand showing some small changes in soil properties. The details of the setup can be found in [1]. We used the DBH measured during 5 complete censuses in order to compute annualised growth in mm/year. The raw data of the experiment was manually rearranged into six files with LibreOffice Calc, each file corresponding to a complete census. We computed annualised growth of each tree in mm/y as the difference between two consecutive censuses divided by the time between the two censuses. In case of mortality of the tree between two censuses, the data was discarded. We computed the neperian logarithm of diameter and of growth (with a constant for growth in order to avoid undefined values).
The final dataset included 64125 DBH measurements corresponding to 13531 trees.
Figure 1.1: The original data without negative growth values (A) and the log-transformed growth (B).
Figure 1.1 shows the distribution of growth after removing negative values, and the latter with log-transformed values.
Figure 1.2: Design of a plot. Each point is a tree and the associated number is the tag of the tree.
Figure 1.2 shows the disposition of the trees in a single plot. There are 14 clones times 10 repetitions, so 140 plots with this same design.
2 Competition index
Figure 2.1: Plot of the growth versus the diameter. Each colour represents a tree age.
Figure 2.1 shows the age of the trees has a big influence on the values of growth but also on the relationship between growth and diameter : the slope is smaller with time, indicating that for the same diameter, growth is slower through time. This is likely an effect of competition for light and possibly underground resources, since as the trees grow their capacity to capture resources increases. Therefore, we computed a competition index to integrate this effect in the growth model. The competition index is computed for each tree which is not on the edge of a plot. It is the sum of the basal areas of the 8 direct neighbours divided by the area of the rectangle that comprises all the direct neighbours. It was then log-transformed.
\(C_{i, t} = \sum BA_{neighbours(i, t)}\)
3 Statistical growth model
In order to partition the variance of individual growth data, we built a hierarchical Bayesian model and used Stan and the package brms [2,3] in Rstudio to implement it in R. We used the following parameters : n.adapt = 1000 ; n.burn = 1000 ; n.iter = 5000 ; n.thin = 5. Our model incorporated a fixed effect on the intercept (\(\beta_0\)), on the slope of diameter D (\(\beta_1\)),and on the competition index C (\(\beta_2\)) and several random effects, namely temporal (date of census, \(b_t\)), individual (tree identifier, \(b_i\)), spatial (block, \(b_b\)), and genotype (\(b_g\)).
Variables are scaled to help the convergence of the model.
\(ln(G_{it}+1) = (\beta_0 + b_i + b_b + b_g + b_t) + \beta_1 \times ln(D_{it}) + \beta_2 \times ln(C_{it}) + \epsilon_{it}\)
Priors
\(\beta_0 \sim \mathcal{N}(mean=0, var = 1), iid\)
\(\beta_1 \sim \mathcal{N}(mean=0, var = 1), iid\)
\(\beta_2 \sim \mathcal{N}(mean=0, var = 1), iid\)
\(b_i \sim \mathcal{N}(mean=0, var=V_i), iid\)
\(b_b \sim \mathcal{N}(mean=0, var=V_b), iid\)
\(b_g \sim \mathcal{N}(mean=0, var=V_g), iid\)
\(b_t \sim \mathcal{N}(mean=0, var=V_t), iid\)
\(\epsilon_{it} \sim \mathcal{N}(mean=0, var=V), iid\)
Hyperpriors
\(V_i \sim \mathcal{IG}(shape = 10^{-3}, rate = 10^{-3}), iid\)
4 Results of the model and variance partitioning
After convergence of the model (checked visually), we examined the variance of each random effect, and this enabled us to perform a variance partitioning.
## No divergences to plot.
Figure 4.1: Trace of the posteriors of the infered parameters
Figure 4.2: Density of the posteriors of the inferred parameters
Figure 4.3: Trace of the temporal random effects
Figure 4.4: Trace of the genotype random effects
Figure 4.5: Trace of the spatial (block) random effects.
Figure 4.6: Mean values and 95% confidence interval of the temporal, genetic and spatial and random effects.
| Intercept (\(\beta_0\)) | Diameter (\(\beta_1\)) | Competition (\(\beta_2\)) | Individual variance (\(V_i\)) | Block variance (\(V_b\)) | Genetic variance (\(V_g\)) | Temporal variance (\(V_t\)) | Residuals variance (\(V\)) | |
|---|---|---|---|---|---|---|---|---|
| Estimate | -2.5e-02 | 5.5e-01 | -2.7e-01 | 2.3e-01 | 5.5e-02 | 1.3e-01 | 1.2e+00 | 5.1e-01 |
| Estimation error | 4.6e-01 | 5e-03 | 9e-03 | 4.1e-03 | 1.6e-02 | 2.7e-02 | 5.5e-01 | 2e-03 |
We found that the two most important contributors to variance were the date and individual identity. High estimation error for the intercept and the temporal random effect must be noted. The proportion of variance represented by each random effect and the residual variance are computed to visualise the variance partition.
Figure 4.7: Proportion of each variance component of the unexplained variance.
The model showed that individual tree growth was a function of tree size and competition with neighbouring trees, and that variance around this model was mostly due to a temporal effect as well as an individual effect (Table 4.1, Figure 1.2). The effect of the genotype was quite small, and the effect of block was even smaller. The temporal random effects declined with the date (Figure 4.6, panel A), showing that the effect of the date on growth is negative (the older the trees become, the less they can grow). We attribute this tendency to competition. Therefore, the competition index C did not fully capture the effect of competition on growth. Another explanation is that the diameter slope was not able to fully capture the decrease of tree growth with size, maybe due to geometrical effects of distributing growth around increasing diameter or physiological constraints linked to height. The temporal effect explained the highest fraction of variance (4.1, 1.2). This is due to the negative effect of competition for light, water, and/or nutrients on growth, which increases with the growth of the trees planted at high densities, and to physiological changes occurring with age. Importantly, variability between individuals is much higher than between genotypes. This shows that there is individual variability even if the trees are clones and, therefore, that individual variability is only caused by exogenous factors. This individual effect can be due to the micro-environment where the tree is positioned, but also to some individual history, such as seedling manipulation and plantation. The genotype explained only 6.1% of the variance. Therefore the impact of individual identity on growth is stronger than the effect of genotype. The block had the littlest impact with 2.5% of the variance explained. This means that the physical environment between blocks is quite homogeneous, or that the physical environment does not play a big role in growth. As the experimental design aimed at minimizing environmental variations and selected productive genotypes able to accommodate several environmental conditions [1], this dataset is a strongly conservative test case for our hypothesis.
5 Conclusion
Overall, we found that there is IV within clonal tree plantations. This shows that the environmental factors (in a broad sense : not genetic) have a strong impact on growth and that IV can indeed emerge within a clone.
6 Code implementation
The whole analysis was conducted using the R language [4] in the Rstudio environment [5]. The tables were made with the kableExtra package [6], the figures with the package ggplot2 [7], and the code uses other packages of the Tidyverse [8] (dplyr [9], lubridate [10], magrittr [11]) and other R packages (here [12], bayesplot [13]). The pdf and html documents were produced thanks to the R packages rmarkdown [16], knitr [19] and bookdown [20].